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Abstract 

We analyze the frustrated Heisenberg antiferromagnet defined on a hon- 
eycomb lattice using a Schwinger-boson mean-field theory. The spin-wave 
velocity and the susceptibility are presented as functions of the strength of 
the frustrating interaction for spin S = ^, and the dynamic structure factor 
is calculated for various temperatures and frustrations. For large 5, we find 
an increased Neel stability with respect to the classical case. 
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INTRODUCTION 



The two-dimensional Heisenberg antiferromagnet has received much attention during 
the last years. One of the reasons for this new interest stems from the discovery of high-T^ 
superconductors, which has encouraged comprehensive investigations of the model defined 
on a square lattice. In particular, there is now convincing evidence that it has a Neel ordered 
ground state for nearest-neighbor interactions even for spin S* = | [|I|. For the honeycomb 
lattice, however, the situation is less clear. The special features of this lattice make it more 
intricate and it has so far been paid less attention. 

The honeycomb lattice is bipartite and has a Neel state as its classical ground state. 
However, since it has the smallest possible two-dimensional coordination number {z = 3), 
quantum fluctuations are expected to be larger than for a square lattice. Hence, it can be 
used for analyzing whether quantum fluctuations may destroy antiferromagnetic long-range 
order (AF LRO) in higher dimensions than one. There are strong indications though that 
this is not the case, even for S = ^. Extrapolations of finite-lattice diagonalizations and 
of Monte Carlo simulations , as well as first-order corrections to spin- wave theory [^] and 
series expansions around the Ising limit all give a finite, albeit much reduced, Neel-order 
parameter. The AF LRO at T = is also found in our approach. 

Another property of the honeycomb lattice is that it is not a Bravais lattice, i.e., trans- 
lational invariance of the full lattice is broken for any type of state. Hence, for a transition 
from magnetic disorder at finite temperatures to Neel order at T = 0, the spatial symmetry 
is not reduced as for the square lattice. Our treatment here should therefore be possible to 
extend to a Kagome lattice, which is neither a Bravais- nor a bipartite lattice. The non- 
Bravais character has also led to a more exotic speculation, namely that the model possesses 
a Haldane gap if it fails to be Neel-ordered 0. However, as just mentioned, several results 
are not in favor of this possibility. 

After our corroboration of the stability of the Neel state on the honeycomb lattice at 
T = 0, it is of course interesting to consider the behavior under frustrating interactions. 
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How will the thermodynamic parameters change, what will happen to the dynamic structure 
factor, and when will the Neel state finally break down? 

In this paper we calculate several thermodynamic parameters and investigate effects of 
frustration of the spin-| case. For this purpose we employ the Schwinger-boson mean- field 
theory (SBMFT), which has proved successful in incorporating quantum fluctuations ^ 
and which two of us have used for an investigation of frustration in the square case In 
particular, we report results for the spin-wave velocity and transverse susceptibility, with and 
without frustration, at zero temperature. We also calculate the dependence of frustration 
for the mass of spin-wave excitations at low temperature. 

While the weak-frustration limit is well-defined and well-behaved both in spin-wave the- 
ory and SBMFT, the case of strong frustration close to the classical critical frustration is 
more complicated. For the square lattice, for example, the location of the Neel boundary 
has been an issue of considerable interest and it is not settled yet. In the present paper 
we present results for the honeycomb Neel boundary using SBMFT, linear spin-wave theory 
(LSWT), and a mapping to the nonlinear a model with SBMFT parameters calculated here 
for 5 = 1. The qualitative differences between the various methods are the same as for the 
square lattice. 

Finally, we present calculations of the dynamic structure factor, which allows for a de- 
tailed study of the nature of the spin-wave excitations. In particular, we probe the exci- 
tations around the antiferromagnetic Bragg peak and their dependence on frustration and 
temperature. 

THE MODEL 

The honeycomb lattice can be viewed as composed of two interlacing triangular sublat- 
tices. Each site has its three nearest neighbors on the other sublattice and its six second- 
nearest neighbors on its own sublattice, whereas the third-nearest neighbors are again on 
the opposite sublattice. With positive exchange couplings, the nearest-neighbor interac- 
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tions classically lead to a Neel-ordered ground state, which is frustrated by the interactions 
within the sublattices. In our study of the frustrated quantum Heisenberg antiferromagnet, 
we consider the J1-J2 model 

(1) 

r,QL R,/3 

with Ji, J2 > 0. The first sum runs over all sites r of one sublattice and over the three 
nearest-neighbor vectors ex.. The next sum, on the other hand, runs over all lattice sites R 
and three of the six second- nearest-neighbor vectors (i. We have chosen 

1 

«i = —ex + -ey , (2a) 

"2 = — + -Qy , (2b) 

0:3 = -ey , (2c) 

and 

Pi = -^e^ - ^^y > (3a) 

/32 = + ^ey , (3b) 

/33 = -V^e^ , (3c) 



and units where the (nearest-neighbor) lattice constant equals 1. These vectors are depicted 
in Fig. 

The Schwinger-boson approach for the honeycomb lattice is similar to the one on the 
square lattice. However, as the former is not a Bravais lattice, complex phases enter geomet- 
ric factors and it is necessary to treat the sublattices separately. In the following, we give 
the main formulas for SBMFT on a honeycomb lattice, which are to be compared with the 
corresponding formulas for the square case (Eqs. (4)-(12) in Ref. [§). The spin operators 
Sfi at each lattice site are replaced by two species of Schwinger bosons foj^jj (/^ = 1,2) via 

Sr = ^^U'^^.Ah > (4) 
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with the local constraints fcj^jj&^j^ = "^S. Here cr — {a^, a^, a^) is the vector of Pauli matrices, 
and summation over repeated (Greek) indices is implied. The spins on one sublattice are 
rotated by tt around the y-axis, which leads us to the following Hamiltonian: 

r,a R,I3 R 

where we have included the local constraints with Lagrange multipliers Xr at each site. The 
expression for the summands is 

^H,s = l-^ks^ny-S' , (6) 
with Xj^ g any of the two link operators 

^R,6 = ^IrKr+S + ^R^R+S ■ (7b) 

The mean-field theory is finally generated by the Hartree-Fock decoupling 

■ '^R,6'^R,6 ■ ^ '^R,s{'^R,s) + {'^R,s)'^R,S " {^R,s) i'^R,d) ' (§) 

where the link fields Qi = (^^^c«) Q2 = {^r^) taken to be uniform and real. In our 
mean- field treatment, we replace the local Lagrange multipliers Xr by a single parameter A. 
After Fourier-transforming the Schwinger bosons independently on each sublattice. 



b 



IxR 



b 



HR' 



with R and R on different sublattices, we use a (complex) Bogoliubov transformation 

%k = cosh ^fcc^fc + sinh , (10a) 

%u = ^""^ cosh i?fec;^ + sinh , (10b) 

to diagonalize the resulting Hamiltonian by choosing and dk properly. This yields free 
bosons with dispersion relation 
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OJk= ^{\ + ^J2Q2l2,kf-\PlQlllA^ , (11) 

where the geometrical factors are given by 

72,fc = ^$^cos(fc-/3) . (12b) 
The mean-field equations are the minimization equations for the free energy: 

j ^ cosh(2^9fc) (nfc + i) - (5 + 1) = , (13a) 
j ^sinh(2i?fe)|7i,fc|(nfc + 1) - ^Qi = , (13b) 
j ^ cosh(2?9fc)72,fc(^fc + I) - = , (13c) 

with A = 87rV(3V3) the area of the reciprocal unit cell of a sublattice, nu = [exp^PuJk) — 1] ^ 
being the Bose occupation number and 

Solving the mean-field equations yields values for Qi, Q2, and A, and also the opportunity 
to determine thermodynamic quantities at finite temperatures. 

In two dimensions there can be LRO only at zero temperature. This means that exactly 
at T = there is an abrupt phase transition, and we must be careful when reducing T to 
zero. As the AF LRO corresponds to a condensation of the Schwinger bosons, we obtain 
the mean-field equations at T = from Eqs. (|l^) by the replacement 

!}h^S*5^^\k) , (15) 

where S* is a new unknown quantity measuring the Bose condensate. On the other hand, 
we can now determine A from Qi and Q2, as there is no gap in the spin-wave spectrum 
{uk=o = in (O)) at T = 0. 
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THERMODYNAMIC PARAMETERS 



For small k the dispersion relation ([TT|) takes a relativistic form: = (mcy + |fcp, 
with the spin-wave velocity c being 

C=l^jJ!Ql-2J2Q2{X + 3J2Q2) , (16) 

and a mass m in the energy gap A = mc^ of the spin-wave excitations: 

A = = y/(A + 3J2Q2)2-(|Jigi)2 . (17) 

We have calculated the spin-wave velocity c at T = with the zero-temperature formal- 
ism and by extrapolation from finite T. The results of the two methods are the same within 
the accuracy of the extrapolation, and they are shown in Fig. |] together with the result for 
LSWT. Without frustration one can compare our result with the first-order SWT (FSWT) 
obtained from l/S-corrections to LSWT. In fact, they give exactly the same result 

c = 1.2098417 * clswt = 1.2832309 Ji . (18) 

Next we consider the gap, or the mass m. As in the square case [jl|, without frustration 
it behaves like m oc exp(— const/T) as T ^ 0. For fixed temperatures the mass (=the gap) 
varies along with frustration. As can be seen in Fig. |^, it increases nearly exponentially with 
frustration, in analogy with the square case. 

Our last thermodynamic parameter is the uniform transverse susceptibility x±- To cal- 
culate this, we first obtain the rotational average x finite T from the static structure 
factor Sxx{q = 0, t = 0) by 

X=-^Sxxiq = 0,t = 0) , (19) 
Kb! 

where 

5.x(q = 0, t = 0) = ^ I ^infc(nfc + 1) . (20) 
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Here, the factor | has been inserted on the right-hand side of ( pOD to make up for the mean- 
field treatment of the local constraints 0. Second, we extrapolate x down to T = 0, and 
use that at T = the longitudinal susceptibility x\\ = 0, so that x± = |x- The result 
is shown in Fig. ^, where one can also see that the LSWT is independent of frustration. 
The extrapolation scheme reduces the precision of the T = values. For a test, we have 
evaluated the exact SBMFT result at T = with no frustration by using a formula obtained 
for the square lattice with integrals calculated for the honeycomb lattice 0]. It reads 

= O.3997332/6J1 , (21) 

to be compared with our extrapolation x± = (0.41±0.01)/6Ji. Two other results to compare 
with are the FSWT resuh [|] x± = 0.274/6Ji and the series resuh @] x± = 0.454/6Ji, which 
shows that the SBMFT gives a result close to the one obtained by a series expansion [Q. 



STABILITY OF THE NEEL STATE 

The classical Neel order vanishes as the frustration a = J2/J1 reaches the critical value 
Oci = 1/6, above which each triangular sublattice becomes ordered and totally decoupled 
from the other sublattice. 

The LSWT correction to the classical Neel boundary for the square lattice was obtained 
by Chandra and Dougot JTOf by finding the frustration at which the sublattice magnetization 
vanishes. A repetition of their calculation reduces the Neel stability also in the honeycomb 
case, as shown in Fig. ^ The square-case reduction was soon challenged by modified spin- 
wave theory calculations (a theory very similar to SBMFT) in Refs. [jll[ and by SBMFT 



results [|T^, which instead pointed to an increase of the Neel region. The discrepancy between 
SWT and SBMFT-like theories has later been argued to exist only in LSWT, and should 



vanish in FSWT through the cancellation of two logarithmic divergences |]T3|-0 . Indeed, for 



large 5, FSWT and SBMFT give the same magnetization [|I^. There seems to be general 
consensus that for large S there is a stabilization of Neel order beyond the classical Neel 
boundary: quantum fluctuations stabilize a state which is classically forbidden. 
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We have also calculated the Neel boundary by using our T = mean-field equations. 
This was done by determining the frustration at which the Bose condensate amplitude 5**, 
and hence the sublattice magnetization, vanishes. This curve is also presented in Fig. ^. As 
in the square case, it is stabilized beyond the classical region, and enhances the Neel stability 
even for S = ^. However, in this latter limit of 5 = |, the phase boundary obtained by the 
large-S* theories can be trusted less. The most extensive exact diagonalizations of the Ji — J2 
model with extrapolations to the thermodynamic limit indeed reduces the Neel stability [ll6| . 

In another approach, we use a mapping to the nonlinear a model with quantum corrected 
thermodynamic parameters for the frustrated spin-i antiferromagnet on a square lattice 



||17|,|8[|. The method considerably reduces the Neel stability in that case. In Ref. IjT^ the 
formulas for the honeycomb lattice were also given, but at that time the quantum corrections 
were not known. Using the values obtained here, these formulas yield that already without 
frustration the bare coupling qq of the nonlinear a model exceeds the fixed point value Qc, 

-go = [SnVsy/'^^^^^^ ^ 13.64 > = An ^ 12.57 , (22) 

CSBMFT XSBMFT 

assuming the hydrodynamic relation p = x(? ■ This means that the model is in its massive, 
disordered state, and that for S = \ there is no Neel order even without frustration. This 
result contradicts the ones obtained previously, and most likely this method underestimates 
the Neel stability. Therefore, the result obtained in Ref. § for the critical frustration on 
the square lattice is most probably too low as well, although we believe the tendency of a 
reduced Neel stability for S* = | to be valid both on the square and the honeycomb lattice. 

To sum up, qualitatively SBMFT and LSWT give the same Neel boundaries as on the 
square lattice, while the effective-action method underestimates the Neel stability for S = 
|. For 5* large, the SBMFT phase boundary probably can be trusted, but much work 
remains to be done to determine the disordering phase transition for the S = ^ honeycomb 
antiferromagnet . 
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THE DYNAMIC STRUCTURE FACTOR 



Finally we examine the spin dynamics by means of the dynamic structure factor S{q, cu). 
For T 7^ there is no LRO, but local Neel order still supports spin- waves with not too- long 
wavelengths. For finite temperatures, the structure factor is rotationally invariant and can 
be written as 

suq, t)^Yl ^''-^''-'''Hsmshm , (23) 

where R and R' run over all lattice sites. Notice that since R — R is not always a multiple 
of the lattice vectors /3, this is not a spatial Fourier transform as it is for a Bravais lattice. 
The (time) Fourier transform is a sum of two terms 

S,x{q, u;)^ Jdt e'^'SUQ, t) = Si{q, u) + S^iq, u) , (24) 
which evaluate to 

Si{q, u) J (-^fe-fe+^g + ^^ rikiuk+q + l)S{u;k - Uk+q + ^) , (25a) 

S2{q.uj) = \j ^ ~ {nu + &{uj)){nk+q + e{uj))5{ujk + uj^+q - W\) , (25b) 
where 

fk,k+q = cosh(2i?fc) cosh(2i?fc+q) - cos(v3fc - (fik+q) sinh(2t?fc) sinh(2i?fc+q) , (26) 

Q{u) is the Heavyside step function, and ip^ is the phase of 71, fc- The non-Bravais character 
of the honeycomb lattice is reflected by the non-zero value of (fk- 

The two terms <Si and S2 correspond to different physical processes and give rise to two 
distinct peaks in the structure factor. In neutron-scattering terms, the first term Si describes 
a simple scattering process of a neutron against a spin- wave, while the other term, ^2, cor- 
responds to a scattering process where two spin waves are created or annihilated. (However, 
this double-peak structure may be difficult to resolve experimentally at low temperatures. 
At higher temperatures the ^2 peak may be too small and smeared out to be detected.) 
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Our results are presented in two series of plots in Figs. ^ and ^ for q close and parallel 
to the Neel ordering vector 



(Compare with plots for the square case in Ref. [|T2|.) All (but one) contain two "peaks", 
where one stems from Si and the other from S2 as indicated in Fig. ^a. 

In Fig. H, T = I Ji (constant) and we vary frustration (||a: J2/ Ji = and Hb: J2/ Ji = 0.1 
< critical frustration). Near the ordering vector we get "overdamped" spin- waves that turn 
into a magnetic Bragg peak as T is lowered. Outside this region, we have the normal spin- 
wave dispersion (5 oc 6{u! — c\q\)). Note that S2 has a small gap in the dispersion. (It 
is actually a double gap, as two spin waves are created.) The gap increases along with 
frustration, as can be seen when comparing Figs, ^a and ^d. (The increase of the gap with 
frustration has also been depicted in Fig. |^.) Above the gap, SWT-type dispersion is back 
to normal. 

In Fig. 1^, J2/J1 = 0.22 > the critical frustration. Being on the disordered side of the 
critical frustration, there is no Bragg peak at T = 0. This remedies the divergencies in the 
structure factor and makes it possible to calculate the structure factor at zero temperature 
within our formalism. In Fig. 0a the results for T = are shown. Here, the iSi-peak has 
vanished since there are no thermally excited spin-waves to scatter against. The ^2 process 
still gives a peak due to creation of spin-wave pairs. At a temperature increase, the Si peak 
reappears as can be seen in Fig. 



CONCLUSIONS 

SBMFT works well for calculating thermodynamic properties for the Heisenberg anti- 
ferromagnet on a honeycomb lattice, although this is not a Bravais lattice. For S = ^, the 
spin-wave velocity and the transverse susceptibility both decrease slowly with frustration, 
while the excitation gap grows nearly exponentially with increasing frustration. For large S 
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AF LRO is found to be stabilized beyond the classical boundary. The case of small S calls 
for further investigation. Finally, the dynamic structure factor shows a double-peak struc- 
ture with a gap for finite temperatures. In general, the qualitative behavior is the same as 
for the square lattice. 
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FIGURES 

FIG. 1. Honeycomb lattice with one of the two triangular sublattices marked. The vectors a 

and /3 point to the nearest neighbors and half the number of second-nearest neighbors of a site, 
respectively. 

FIG. 2. The T = spin- wave velocity vs frustration J2/J1 for S* = i. The full line shows the 
SBMFT result and the dashed line the LSWT result. 

FIG. 3. The logarithm of the mass vs frustration for SBMFT at temperature T = ^Ji and 

s = l 

FIG. 4. The T = uniform transverse susceptibility vs frustration for S = ^. The full and 
dashed lines correspond to SBMFT and LSWT, respectively. 

FIG. 5. The inverse of the critical spin Sc as a function of frustration. The lines show the 
LSWT and the SBMFT result. They coincide for J2 = and for J2/J1 = 1/6 at the classical point 
S = 00. Note the Neel stability beyond the classical value in the SBMFT case. 

FIG. 6. The dynamic structure factor at a low temperature T = iji close to the ordering 
vector qq. In (a) the frustration is J2/J1 = 0, and in (b) it is J2/J1 = 0.1 < critical frustration. 
The peak corresponding to S2 has a gap (starts at finite u) which increases with frustration. 

FIG. 7. The dynamic structure factor for J2/J1 = 0.22, which is a value above the critical 
frustration. In (a) the result for T = is shown, while in (b) we have T = ^Ji- Note that the Si 
peak has vanished in (a) due to the lack of thermally excited spin waves. 
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